Electronic states near a quantum fluctuating point vortex in a d-wave superconductor: 

Dirac fermion theory 
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We introduce a simple model of the low energy electronic states in the vicinity of a vortex under- 
going quantum zero-point motion in a d-wave superconductor. The vortex is treated as a point flux 
tube, carrying n flux of an auxiliary U(l) gauge field, which executes simple harmonic motion in a 
pinning potential. The nodal Bogoliubov quasiparticles are represented by Dirac fermions with unit 
U(l) gauge charge. The energy dependence of the local density of electronic states (LDOS) at the 
vortex center has no zero bias peak; instead, small satellite features appear, driven by transitions 
between different vortex eigenmodes. These results are qualitatively consistent with scanning tun- 
neling microscopy measurements of the sub-gap LDOS in cuprate superconductors. Furthermore, as 
argued in L. Balents et al, Phys. Rev. B 71, 144508 (2005), the zero-point vortex motion also leads 
naturally to the observed periodic modulations in the spatial dependence of the sub-gap LDOS. 



I. INTRODUCTION 

In a recent paper-, hereafter referred to as I, two of the 
present authors considered influence of the vortex zero- 
point motion on the energy dependence of local density 
of electronic states (LDOS) in a s-wave superconductor. 
Here we will present a study of vortex zero point motion 
in two-dimensional d-wave superconductors, with an eye 
to application to the cuprate supercondcutors. A direct 
application of the methods of paper I to d-wave super- 
conductors is discussed in AppendixEl but this approach 
has significant limitations. It makes a gradient expansion 
in the spatial dependence of the gap function and so can- 
not be applied in the limit of small coherence length, 
and is not designed to efficiently extract the important 
effects of the low-energy nodal quasiparticles. The body 
of the present paper will describe a new approach which 
overcomes these limitations. 

Our focus on the zero-point motion of vortices in the 
cuprates is motivated by a previous proposal that peri- 
odic modulations in the spatial dependence of the LDOS 
inevitably appear over the region the vortex executes 
its quantum zero-point motionSi^. Such periodic LDOS 
modulations have been observed in scanning tunneling 
microscopy (STM) studies of the vortex in the cuprate 
superconductors^. Here we study whether the vortex 
motion could also explain the energy dependence of the 
LDOS near the vortex core. 

A solution of the Bogoliubov-de Gennes (BdG) equa- 
tions for a vortex in a d-wave superconductor leads to 
a large peak, as a function of energy, at zero bias at 
the vortex cente r^i^i^ . No such peak is observed in the 
experiments. The inclusion of additional density wave 
ordergSiiS. can suppress the zero bias peak in favor of 
spectral weight at satellites, but only under conditions in 
which the strength of the order parameter is unaccept- 
ably large. The theories require nearly complete mag- 
netic/charge ordering, and there is no indication from 
e.g. neutron scattering experiments that such a large or- 
dering can be present in the optimally doped Bi-based 



cuprates. 

Paper I extended the BdG equations to include vortex 
zero-point motion for s-wave superconductors: it found 
that zero bias peak was reduced (but not eliminated), 
with a transfer of spectral weight to energies of order 
the vortex oscillation frequency. This was done in an 
approach that performed a gradient expansion in the gap 
function, which is only valid for a large core size. 

The present paper examines an alternative approach to 
computing the sub-gap energy dependence of the LDOS 
near the vortex. Our starting point is an observation 
by Tsuchiura et ali^ that the zero bias peak disappears 
when the vortex core size becomes of order the inverse 
Fermi wavevector. Consequently, we will examine a con- 
tinuum model in which the vortex core becomes point- 
like, with a vanishing radius. In this situation, no gradi- 
ent expansion can be performed, and it is essential to 
compute the effects of vortex zero-point motion with- 
out expanding in powers of the vortex position. We 
will demonstrate, within the context of a simple model, 
that not only is the zero-bias peak eliminated, but the 
vortex motion leads to satellite features associated with 
transitions between different vortex vibrational states. 
We will argue that these features are appealing candi- 
dates for explaining the sub-gap peaks observed in STM 
experiments^i^iiSiiii^. 

In a sense, the analysis of this paper can be viewed 
as a method of computing the influence of "phase 
fluctuations'— on the LDOS. However, instead of ex- 
plicitly integrating over the phase degrees of freedom, we 
represent them by a collective co-ordinate, the position 
of the vortex. Apart from the benefit of the compact rep- 
resentation, this allows us to physically interpret the cou- 
pling of the phase fluctuations to the electronic quasipar- 
ticles. Furthermore, inclusion of Aharanov-Bohm phase 
factors exposes the subtle relationship of the vortex fluc- 
tuations to competing order parameters^. 

The outline of the paper is as follows. Our model will 
be introduced in Section m along with an initial discus- 
sion of its characteristic properties. The details of our 
calculations are presented in Section IIIII After setting 
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up the formalism, we first analyze the LDOS influenced 
only by the vortex zero-point quantum motion, and then 
include the resonant scattering at the lowest order of per- 
turbation theory. Plots of the full LDOS are shown in the 
Sections IIII CI and IIII Dl Wc summarize our results and 
their relation to experiments in Section llVI 



II. THE MODEL 

We consider a single vortex coupled to quasiparticles 
in a clean d-wave superconductor at zero temperature. 
The vortex is assumed to experience a harmonic trapping 
potential in which it can carry out its quantum zero-point 
motion. Such a potential can result from interactions 
between vortices in a vortex lattice, or from a pinning 
impurity. The Hamiltonian can be written as {h = 1): 

+ J2 J d\¥ir)HBdG{r)^ (r) . (1) 

nodes 

For now we neglect the Magnus force on the vortex, as- 
suming that it is much smaller that the trapping force; 
we shall consider it carefully in the Section llllDI The op- 
erators and Tv are the canonical momentum and posi- 
tion of the vortex. For simplicity, we place the vortex in 
an isotropic trap characterized by a single harmonic fre- 
quency LUv It has been shown that presence of static vor- 
tices does not qualitatively change the low-energy spec- 
trum of nodal d-wave quasiparticlesSiiS^, and in this 
paper we will find that the gapless nodes also survive 
quantum fluctuations of vortices. Therefore, quasiparti- 
cles are massless Dirac spinors and we describe them by 
Nambu operators '^^{r) and ^'(r) defined at every node 
in the c?-wave spectrum. Being interested only in the low- 
energy dynamics, we linearize the Bogoliubov-de Gennes 
Hamiltonian in the vicinity of gap-nodes, and apply a 
convenient Franz- Tesanovic unitary transformation. For 
the node p « kfX we haveiSiii: 



BdG 



VfiPx + ax)(T^ + t^A(Pa + o-v)^^ + mvfVxI , (2) 



and we need not explicitly worry about the other nodes 
because their linearized Hamiltonians are related by uni- 
tary transformations. Note that we formally work with 
dxy symmetry, which is related to d^i^yi by a rotation. 
In this expression, p is the quasiparticle momentum op- 
erator relative to the node, m is the electron mass, Vi and 
v/^ are the Fermi and gap velocities respectively, and 
are Pauli matrices. The effective gauge field a is pro- 
portional to the phase gradient of the superconducting 
order parameter, and thus corresponds to the 7r-flux cen- 
tered at the vortex location. The role of a is to im- 
plement the statistical interaction between vortices and 
quasiparticles, so that when a quasiparticle completes a 
circle around the vortex its wavefunction changes sign. 



The supercurrent velocity field v appears in the Hamil- 
tonian because the supercurrents that circulate around 
the vortex give rise to Doppler shifts of quasiparticle en- 
ergies. Note, however, that only the projection of v on 
the nodal direction matters, and also that Doppler shift 
effects decrease rapidly beyond the London penetration 
depth from the vortex center. 

Our goal is to elucidate only the qualitative features 
of quasiparticle spectra in the vicinity of a fiuctuat- 
ing vortex. Hence, we will make several simplifications 
that are not quantitatively justified in realistic situations, 
but the obtained spectra will nevertheless be remark- 
ably similar to the spectra observed in experiments. The 
complexity of calculations is considerably reduced if the 
Doppler shift and all effects of anisotropy are neglected. 
The Bogoliubov-de Gennes Hamiltonian near the node 
p K, ifefi becomes: 



H 



BdG 



Px 



Py 



Py + ay -{px + ax) 



(3) 



where we have set Vf and wa to unity. Wc will work in 
the units h = V{ — 1 at all times, except when we discuss 
scales. 

The coupling between the fermions and the vortex 
arises from the gauge field. This is specified by the in- 
stantaneous position operator of the vortex: 



a(r,rv) 



z X (r - Tv) 
2|r-rvP 



(4) 



assuming that we can regard the vortex core to be negli- 
gibly small. 

The physical reason for these simplifications comes 
from our expectation that it is the vortex quantum fiuc- 
tuations that produce the sub-gap peaks at about 7 
meV from the Fermi level in the quasiparticle density 
of statesiSiiiii. Even if we regard the core of a static 
vortex to be infinitely small, as was done in the expres- 
sion the zero-point quantum motion of the vortex 
in the harmonic trap creates a finite central region in 
which the supercurrents are gradually suppressed. The 
circulating superfiow is then the strongest at some finite 
radius. Quasiparticles can scatter from this supercurrent 
distribution through any number of virtual excited states 
of the oscillating vortex. The result is a weak resonant 
state, a metastable binding between the vortex and a 
quasiparticle. In our simplified calculation this effect is 
solely due to the electron Berry phase, but the Doppler 
shift should also contribute it, even more directly as a 
partially attractive potential. Indeed, our calculation will 
qualitatively reproduce the sub-gap peaks in the density 
of states, which are reminiscent of a localized metastable 
state in the vortex core, and in our model originate from 
resonant scattering. In the light of this scenario, we do 
not expect anisotropics to play a qualitatively significant 
role. The Doppler shift may be important in enhanc- 
ing these effects, but the mere demonstration of the sub- 
gap peaks through the Berry phase alone makes the case 
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for relevance of the vortex quantum motion. Another 
possible mechanism is resonant scattering from the con- 
ventional vortex core, where the superconducting order 
parameter magnitude is suppressed. However, in rf-wave 
superconductors the core is very small (a few lattice spac- 
ings), so that the corresponding resonance would occur 
only at relatively large energy. Note that we are not con- 
sidering the possibility of bound core states: all quasipar- 
ticle states are expected to be extended in pure d-wave 
superconductors, even in the presence of vorticesSiiSiii. 

Our main result is the quasiparticle local density of 
states (LDOS) in the vicinity of a quantum oscillating 
vortex. We calculate the LDOS perturbatively, using 

a^(^y (5) 

as a small dimensionless parameter. We treat the vortex 
effective mass as an independent quantity, although 
it can be determined microscopicalljii^. After all simpli- 
fications, we can write the LDOS as the following scaling 
function of energy e and distance r from the origin (center 
of the vortex trap): 

*-) = ||E«-^.(£:.|^:«)- w 

I n— ^ ^ 

The universal functions Fn{x,y\a) can be regarded as 
being of the same order of magnitude for any given finite 
arguments. Note, however, that this is not a Taylor ex- 
pansion. The extra dependence of Fn on a is a somewhat 
unusual consequence of the proper choice of perturbation. 
It is usual to start from a non-interacting theory and in- 
clude all interactions as perturbations. Unfortunately, 
interactions between quasiparticles and vortices are not 
weak, and cannot be characterized by a small expansion 
parameter. Instead, it is better to include the influence 
of the vortex zero-point quantum motion on quasiparti- 
cles at the zeroth order, and treat the resonant scattering 
of quasiparticles from the oscillating vortex as a pertur- 
bation. Thus we perform a partial resummation of the a 
expansion, and this accounts for the a dependence of the 
co-efficients in Eq. ^ ; features of the zero-point vortex 
motion, also parameterized by a, appear in the LDOS at 
all orders of perturbation theory. The small parameter 
a is the ratio of such perturbation's energy scale and the 
energy barrier hw^ for virtual transitions. 

As a matter of fact, only the zero-point vortex oscilla- 
tions set the scale for the zeroth order term in ® , which 
results in an additional scaling property of Fo(x, y; a): 

Fo{x,y;a) ^ aFo{x/a,y;l) . (7) 

One should not be misled to conclude that Fq becomes 
small when a — > 0: as we will show, Fq{x, y; a) x/(2Tr) 
for X ^ a. 



III. PERTURBATION THEORY 

The Hamiltonian can be expressed in a fashion 
analogous to the Holstein-Primakoff expansion for the 
spin systems. In absence of quasiparticles the vortex is 
modeled by a two-dimensional linear harmonic oscilla- 
tor, whose eigenstates \nx,ny) are characterized by two 
integer quantum numbers, Ux and Uy. Let us define the 
following "matrix elements" : 

V(„i^,„i„),(„2,,„2„)(^) = {nix,niy\HBdG{r)\n2x,n2y) . 

(8) 

We also introduce creation operators, 6j, and 6^, which 
raise the vortex quantum number: 

bl\nx,ny) = y/ux + l\nx + l,ny), (9) 

and similarly for b^. If we insert the identity operators on 
the left and right side of HbcIG in the equation , and 
resolve them in terms of the vortex eigenstates \nx,ny)^ 
we can systematically write: 

■^^HBdG-^ =*^{l^(0,0),(0,0)+ (10) 

V(i,o), (0,0)4 + h.c. -f V(o,i),(o,o)&J + h.c. + 

(V(i,o),(i,o) - V(o,o),(o,o)) 4^ + 

(V(o,i),(o,i) - V(o,o),(o,o)) blby + 
V(i,o),(o,i)4^y + h.c. + V(^i,i)^(Q,Q)blbl + h.c. + 

^V(2,o),(o,o)&i4 + h.c. + 
^V(o^2)Xofi)hlhl + h.c. + ---]^ . 

It can be easily seen that the quasiparticle momen- 
tum operators p appear only in the lowest order term 
^(0,0), (0,0)- In all clastic processes (such as b^bx) the mo- 
mentum operators that appear in the matrix elements V 
are cancelled out. On the other hand, in all inelastic pro- 
cesses the momentum operators do not appear because 
they do not excite vortex states. Therefore, for practi- 
cal purposes we can replace all matrix elements except 
V(o,o),(o,o) with: 

l'("i.,"i„),(n2x,n2jl'^J - a^(^) -a'^{r) J ' 
where 

a'(r) EE a(„^^,„^j^(„^^^„^j(r) (12) 
= {nix,niy\a{r,r^)\n2x,n2y) . 

Note that the bare effective gauge field a(r,rv), given 
by 10} , does not define any energy or length scale. The 
scales are introduced through the states \nx,ny) of the 
quantum harmonic oscillator. The harmonic oscillator 
wavefunctions depend only on the dimensionless coordi- 
nates ^ = y/rru/^r, so that the energy scale associated 



4 



^) is {hujv'm^vfy^'^ (when all physi- 



with V(„^_^,„^j,(„,^,„,_ 
cal constants are put in their places). 

The expansion H1Q|I separates various processes in 
which quasiparticles can scatter elastically or inelastically 
from a vortex. At low enough energies the expansion can 
be truncated after a few lowest order terms. In the fol- 
lowing we will set up a perturbation theory, with the 
unperturbed Hamiltonian: 



Ho - {blb^ + blby) + 



and the perturbation: 



(0.0), (0,0) 



(13) 



Fi = / dV^-t (1^(1,0), (0,0)4 + l^(o,i),(o,o)fe]; + h.c.) ^ 



(14) 

The more complicated scattering processes in (|1(JII be- 
come qualitatively important only at energies e ~ 2ujv 
and above. 

Before we begin calculations, it is useful to identify 
the small parameter of the perturbation theory. The 
available parameters in the model define two energy 
scales: the energy of vortex quantum oscillations fkOv, 
and the characteristic kinetic energy m^v^. The unper- 
turbed Hamiltonian H13|) is properly characterized by 
the scale \Hq\ = huj^, because the vortex harmonic fre- 
quency sets the energy barrier that needs to be crossed 
in perturbative virtual transitions. The energy scale of 
the perturbation ifT^ is |iJi| = {huj^m^vfY^'^ , as dis- 
cussed above. Therefore, we can define the small param- 
eter a = \Hi\/\Ho\ given by expression (0). The lowest 
order correction pi to the quasiparticle LDOS comes from 
the Fock exchange process (self-energy) , so that pi (x a^. 
Namely, every vertex cx \Hi \ in Feynman diagrams intro- 
duces a factor of a, which becomes apparent after rescal- 
ing all energies by | Hq \ = fiw^ as was done in © . 

In practice, one also has to keep in mind the quasipar- 
ticle cut-off energy A. This third energy scale physically 
comes from the superconducting gap amplitude, beyond 
which linearization of the Bogoliubov-de Gennes Hamil- 
tonian breaks down. Furthermore, renormalization of the 
effective vortex mass due to quasiparticles^® is such that 
m^vf ^ A. Irrespective of its physical origin, numerical 
analysis requires that A be introduced. However, in the 
ideal linearized model that we are concerned with, we can 
in principle take the A ^ oo limit. As long as we treat 
the vortex mass and frequency Wv as independent 
parameters, we find that the spectrum does not depend 
on the numerically introduced A at energies sufhciently 
below A. The functions F„ in © are universal. 



The LDOS due to vortex zero-point quantum 
motion 



We begin by considering the lowest order term in the 
expansion ((BJ for the LDOS: 



(15) 



This term contains effects of the vortex zero-point quan- 
tum motion on quasiparticle spectra, but excludes effects 
of the resonant scattering. For this purpose we will nu- 
merically diagonalize the Hamiltonian H13|) in infinite 
space. 

It will be useful to first understand the solution for a 
static vortex, which is obtained in the limit of nivtOv — *■ 
oo. In this limit there are no vortex fluctuations, so that 
^0,0), (0,0) reduces to (O, with the gauge field given by 
0, and Tv = 0. Exact analytical solutions (properly 
normalized in the infinite space) are then found to bei^: 



-iq- J_;_i (fcr)e*'"^ 
(^^)e*(/-i)0 \ 



(16) 



J, 



iq- Ji+iikr)e 



I < 
I > 



in the representation that diagonalizes the angular mo- 
mentum. Quasiparticles carry the following quantum 
numbers: "charge" q = ±1 (distinguishes particle- like 
and hole- like states) , angular momentum I G Z and radial 
wavevector A: > 0. There are only extended quasiparticle 
states with energies e = q\k\. Ji{kr) are Bessel functions 
of the first kind. The wavefunction for Z = is a linear 
combination of the two forms of ^/j^q in (|16|l . which can 
be uniquely specified only with a separate boundary con- 
dition at the vortex location (an artifact of m^oj^ — > oo 
and a point-like core)ii. The corresponding LDOS poo 
diverges as 1 /r at any energy, due to the I — wavefunc- 
tions, which include J_i{kr): 



cos(2 e r) 


1 1 °° 






/=0 




\e\r < 1 


1 ^ 


\e\r > 1 



.(k|r) 



(17) 



For finite iVy and rriy the vortex fluctuates in a finite re- 
gion whose radius is ~ (wv'tiv)"^''^. This is a property of 
the harmonic oscillator ground-state. The exact quasi- 
particle eigenstates are still characterized by the same 
quantum numbers, and only the dependence on radius r 
is changed. This is so because in the ground-state of our 
model the vortex wavefunction is rotationally symmetric, 
and well localized in a finite region of space. The LDOS 
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po affected only by vortex zero-point oscillations will co- 
incide with H17|l at energies e ^ (mvWv)^^^ or distances 

The first step in diagonalizing is to calculate 



V, 



(0,0), (0,0) 



. From the definition 



we see that Vi 



(0,0), (0,0) 



is given by the expression (jSJ with an effective "gauge 
field" : 



0(0,0), (0,0) (J^) = (0,0|a(r,rv)|0,0) 



(18) 



where a(r,rv) is given by Substituting here the 

ground-state wavefunctions of the linear harmonic oscil- 
lator, we get: 



a(o,o),(o,o)(r-) = (l-e-™^'^"'-') a{r,0) 



(19) 



Details of this derivation are presented in Appendix ^ 
The effective gauge field is not singular at the origin, 
even though the static vortex has a point-like core. This 
indicates that supercurrents are suppressed near the ori- 
gin. Effectively, a "vortex core" has been created, and 
hence the exact wavefunctions are significantly modi- 
fied from H16|l . Behavior of the exact wavefunctions 
-0 = (e**^'~"'^-''^M(r), e''"^u(r)) at small r becomes: 



u{r) cx J;_i(K;r) 



(20) 



and similar for v{r), where ki = ^/k"^ — Imyio^. These 
wavefunctions are always finite or zero at the origin, so 
that the LDOS is not divergent. For r ~> (yj^m^Y^I"^ the 
exact wavefunctions reduce to the form given by (|16|l . 
except that there is a phase shift: tp{r,(j)) = tp°°[r + 
6i^k,4>)- The known behavior for small and large r allows 
one to easily find and normalize exact wavefunctions in 
infinite space: the radial Schrodinger equation is first 
numerically solved in a finite region of space r < tq, and 
then the solution is matched by value and derivative to 
the phase shifted form (|16l) for r > ro ^ (mvWv)"^^^. 

In the figure Qwe plot the LDOS po(e) at the origin 
for several values of the dimensionless parameter . All 
plots collapse to a single curve (shown in the inset) if both 
the energy and LDOS are rescaled by (mvWv)^^^, because 
only this energy appears in the quasiparticle part of the 
Hamiltonian (|13|l . Hence, the following scaling applies: 



-F' 



which in comparison with (|15|l implies: 

Fo{x,y;a) = aFo(x/a,y) = aFQ{x/a,y; 1) 



^ , (21) 



(22) 



At large energies (e ^ (mvWv)^^^), the LDOS at the 
origin begins to reflect the bulk density of states, and 
hence becomes a linear function of energy. Furthermore, 
Po{£, r) in the bulk (r ^ (tovWv)"^^^) must coincide with 
Pooi^if) of the static vortex, so that the expression H17|l 
determines both the large energy and large distance be- 
havior of po(ej r)'- 



X 

2^ 



V?/> 1 



(23) 



Po/(o 




0.5 1 1.5 



e/(Ov 



FIG. 1: The quasiparticle LDOS influenced solely by the vor- 
tex zero-point quantum motion. The LDOS is plotted at the 
origin as a function of energy and o? . All plots can be rescaled 
to a single curve, shown in the inset. 




FIG. 2: The quasiparticle LDOS influenced solely by the vor- 
tex zero-point quantum motion, plotted at the energy e = tJv 
as a function of distance r from the origin. The solid and 
dashed lines correspond to an oscillating {o? = 5), and a 
static vortex respectively. The apparent drop of LDOS for 
the oscillating vortex beyond r ~ 4 is an artifact of numerical 
calculation, which took into account only a small number of 
angular momentum channels (—5 < i < 5). 



Spatial variations of the LDOS are illustrated in the fig- 
ure 121 Zero-point quantum motion of the vortex removes 
the 1/r divergence of the LDOS found in the core of a 
static vortex, while at large distances quasiparticles do 
not feel the vortex fiuctuations. 

Smallness of the vortex core, measured by the coher- 
ence length ^ ^ (tovWv)"^''^, seems to be responsible 
for absence of the zero-energy peak in the quasiparti- 
cle LDOS, as was suggested in the Refini This can be 
inferred by comparing the LDOS obtained here with the 
LDOS of a similar calculation from the AppendixlUl The 
approach from the Appendix works in the opposite limit 
^ > (tovWv)^^^^ and produces a broad zero-energy peak 
in the quasiparticle LDOS. 



B. One-loop correction 

In order to explore the resonant scattering of quasi- 
particles from the fiuctuating vortex, we calculate the 
one-loop correction to the quasiparticle Green's function. 
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This will directly lead to the lowest order LDOS correc- 
tion: 



pi[e,r) = a 



Fi 



er 



(24) 



which we will discuss in Section UlI CI 

The Hamiltonian given by and H14|l defines a 

problem of massless Dirac fermions coupled to a single 
bosonic oscillator. The unperturbed fermion states are 
characterized by "charge" q, radial wavevector k and an- 
gular momentum I. Perturbation theory can be visu- 
alized using the standard Feynman diagram technique. 
The bare propagators (Gq for quasiparticles and Dq'^ for 
the vortex) and vertices are: 



q.kj 



Go{q,k,l,uj) = ^_^(^L,o+) 

^0 \^) — a;2-a;2+iO+ "A'I' 




qi,ll,ki;q2,h,k2 



(25) 

Here, the indices /x and v denote spatial directions x and 
y, and summation over repeated indices will be assumed 
from now on. The vertex operator V^^ ^.^.^^ ^.^ is the 
Fourier transform of V(i o),(o.o) = ^(o.o),(i.o) for /z = a; and 
V(o,i),(o,o) = V(o^o),(o,i) for fJ. = y- Frequency is conserved 
at vertices, but k and I are not. Physical conservation 
of angular momentum is reflected by requirement that 
V'f' , I. .„ , u oc i; the first excited states of the 

2D harmonic oscillator carry angular momentum Z = ±1, 
and in these inelastic scattering processes the quasipar- 
ticle angular momentum changes by one. 

The one- loop self-energy E has two contributions, in 
principle: 



CO -co' 



qi.kiJi 
q2.k2.i2 

Hartree 



^I'Kh q',k',l';co' ^ih'h 

Fock 



(26) 

The Hartree self-energy is zero within the present ap- 
proximations because the vertex that touches the fermion 
loop has the same incoming and outgoing angular mo- 
mentum. The Fock self-energy is: 

dk'—Go{q'J',k',uj') 



X Dl^'^iuj - w')^(,'^,ii,fci;q',i',fe'K','',fc';92,'2,fe2 
_ / ^1^/ qi,h:ki;q' ,1' ,k' ' q' ,1' ,k';q2,l2,k2 

q',1'-' 



LJ - g'(wv + k' - 



(27) 



The one-loop correction to the quasiparticle Green's 
function is given by the simplest diagram that contains 



only one self-energy part: 

G':LM;q.MM^^) - (28) 
Ga{qi,li,ki,uj)Go{q2, h, k2,uj)T.q^j^^kuq2M,k2 i^) ■ 

The remaining task is to calculate the vertex operators 
^qi h ki q2 12 k2 momcutum space. In the original posi- 
tion representation these operators have the form 
where the first-order effective gauge fields (/i e {a;, y}) 



are: 



a" = (l,0|a(r,rv)|0,0) =- 



(29) 



1 2^,2: X I 



(0,l|a(r,rv)|0,0) 



1 2Lzx^ 



e 



Here, ^ = [£,xiS,y) are the usual dimensionless coor- 
dinates of the quantum harmonic oscillator, given by 
^ = ^JrfQJJ^r. For switching to the momentum repre- 
sentation it will be convenient to keep the spinor struc- 
ture and organize various quantities into spinor matrices. 
The role of spin s is taken over by the quantum number 
q in the momentum representation, which can take only 
two values being the sign of energy. First, we define the 
Fourier weight matrix, which is used to translate between 
the position and momentum representations: 



Ti,k{r) 



W+,i,fe(r) U_,;,fc(r) 



(30) 



This matrix is expressed in terms of the eigenfunctions of 
the unperturbed Hamiltonian , which in the position 
representation look like: 



i^q,kA^) 



Uq.Lkir) 
Vq.lA'^) 



(31) 



The quasiparticle spinor operators in position and mo- 
mentum representations are related by: 



*(r-) =^ / dfcTi,fe(r)*i,fc 



(32) 



In the momentum representation, the bare quasiparticle 
propagator is of course diagonal: 



Go{kJ,u) 



Lj-k+iQ+ 






1 

uj+k-ia+ 



(33) 



The momentum representation of the vertex operators 
is obtained in the following manner: 

KM;i2.k2 = / ^'rTlkS^)V'^m^Mi^) , (34) 
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where V^(r) is either V(i_o),{o^) or V(o,i),(o,o) depend- 
ing on ^. When ((TT|l and (|2^ are substituted in this 
expression, it is possible to exactly integrate out the az- 
imuthal angle (p. This implements the physical angular 
momentum conservation, so that V,^ i ; i, = unless 

1^1-^21 = 1. 



C. The LDOS due to resonant scattering 

The local density of states is naturally obtained from 
the quasiparticle Green's function. Including the quasi- 
particles at given position and energy of either spin and 
angular momentum gives the following expression in the 
spinor momentum representation: 



p(e,r) = -isign(e) • Imi ^ J dkidk2 

\ h ,1-2 



(35) 



Tr 



In the absence of perturbations, we would substitute 
here the bare Green's function Gq from It is easy 

to show that the zeroth order term in the LDOS is: 

+ (|z._,,l,|(r)p + |«_,^|,|(r)p)e(-e)} . (36) 

This is the expected result, which could have been im- 
mediately written from the known eigenfunctions H31|l : 
it is plotted as a function of energy in the figure ^ 

The one-loop correction pi(e,r) to LDOS is found by 
substituting (|28|l into H35|l . The correction to the 
Green's function is written in the spinor representation 
as: 

'^/i^fci;/2,fc2('^) = Go(Zi,A:i,w)S,,,fci;;2,fc2(w)Go(?2,fc2,w) , 

(37) 



where the self-energy matrix is: 



(38) 




1 



V; 



Details of the expression for LDOS and numerical proce- 
dures are elaborated in Appendix ^ 

The one-loop correction to LDOS at the vortex center 
is plotted in the figure |21 as a function of energy. Some 
general features can be immediately noted. There is a 
discontinuity of LDOS when the quasiparticle energy is 
equal to the vortex harmonic frequency, e = uJv Also, 
for all values of the small parameter a^, the LDOS has 
a local minimum at zero energy and grows into a peak 
at some finite energy. This feature could be interpreted 
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FIG. 3: The one-loop correction pi to LDOS at the vortex 
center as a function of energy. These plots show evolution of 
pi(e) as the small parameter changes. The magnitude of 
pi scales as a^. Quasiparticle cut-off energy is A = 4ljv 



as formation of a weak metastable state inside the spa- 
tial region spanned by the vortex quantum oscillations. 
Note that the LDOS peak shifts toward lower energies 
as oc rriy decreases: for smaller vortex mass the 
vortex oscillates in a larger region and hence tries to lo- 
calize quasiparticles at lower energies. For very small val- 
ues of some additional features develop in the LDOS, 
such as secondary peaks and dips. It is possible that the 
values oc l/n^ for integer n are special and control 
appearance of various LDOS features (the plots hint to 
this possibility, but are not entirely conclusive). Physi- 
cally, a is the ratio of two length scales: the quasiparticle 
wave-length at the energy e = lu^, and the spatial extent 
of the vortex zero-point quantum oscillations (effective 
core radius). Roughly speaking, features of the LDOS 
change qualitatively when the effective core region grows 
to enclose an additional quasiparticle wavelength l/o^v 

In the figure 2|we plot energy ec of the LDOS sub-gap 
peak as a function of a^. The peak position scales as 
energy of the perturbation that creates it: Ec oc aujv In 
order to see how this result relates to experimental ob- 
servations, we will explore a somewhat different scenario. 
Suppose that the plot in the figure 0]were linear, so that 
ec/u>v were proportional to a^. This would mean that the 
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FIG. 4: (color online) Energy of the sub-gap peak in the 
LDOS as a function of the small parameter o? , and the fit to 
ec oc a (red dashed), ec is defined as the energy of the first 
local maximum in pi(e) from the figure 



sub-gap peak energy ec were proportional to mvvf (see 
equation Q). It has been arguedi^ that the nodal quasi- 
particles renormalize the vortex mass to a value that 
is proportional to the superconducting gap amplitude Aq. 
This renormalization due to quasiparticles is, likely, the 
largest contribution to the effective vortex mass, since 
the usual logarithmic infra-red divergence of the hydro- 
dynamic vortex mass is cut-off by the Coulomb screen- 
ing. Therefore, according to this scenario, Ec oc Aq. STM 
measurements are consistent with such linear scaling of 
the sub-gap peak energy with the superconducting gap 
amplitude. Most notably, the sub-gap peak position was 
experimentally foundi2iiii4 not to depend on the mag- 
netic field, and the considered scenario would agree with 
this observation (even though we work with uj^, which 
could depend on the magnetic field through the inter- 
vortex separation in a vortex lattice). Discrepency be- 
tween our actual result (ec oc hm^Lo^vj) and experi- 
mental observations could be due to a combination of 
large experimental error margins and simplifications in 
our model. Further study is needed to appreciate effects 
of the Doppler shift and other important factors. 

The full quasiparticle LDOS p = po+Pi as a function of 
energy, measured at different distances from the vortex 
center, is plotted in the figure |S1 The sub-gap peak 
height is proportional to inside the core, and gradually 
decreases over a length-scale comparable to the extent of 
vortex zero-point oscillations (mv^v)""'^^^- Below « 
0.2 the sub-gap peak becomes too small to be visible, 
just like all other features of the one-loop correction to 
LDOS. 

Spatial variations in the one-loop correction to LDOS 
Pi(r) are demonstrated in the figure l^jfor two values of 
a^. In both cases pi{r) was calculated at the energy of 
the dominant peak in pi{e;r — 0), which can be seen in 
the figure 13 Correction to the LDOS is largest inside the 
region of the vortex quantum oscillations, and hence in 
agreement with the conclusion that vortex fluctuations 
create a local resonance in the quasiparticle spectrum. 
All spatial variations are symmetric under rotations in 




FIG. 5: Energy scans of the full LDOS p — po+pi at gradually 
increasing distances r from the vortex core. The plots are 
offset vertically for clarity, starting from r = at the bottom, 
and moving up with increments Ar = Q.2uj^^ for — 1 and 
Ar = O.SSoj-^ for = 0.3 (Arwv ^ (Sa)"^). 



our model. 

Finally, in the figure we explore the influence of a 
finite quasiparticle cut-off energy A on the spectra. The 
full LDOS is plotted as function of energy for three values 
of a^. The main consequence of finite A is to reduce the 
density of states at large energies. This effect is more 
pronounced for large when validity of perturbation 
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a =1 , e=Wv 




FIG. 6: The one-loop correction pi of LDOS as a function of 
radius. The plots are taken at energies of the dominant peak 
in /5i(e). Accuracy is lost to the right of the dashed line due 
to the small number of included angular momentum channels. 



theory becomes questionable. In general, the low energy 
LDOS is only weakly affected by A, and the limit A ^ oo 
is well defined. 



D. Influence of the Magnus force on the LDOS 

When a vortex moves with respect to the superfluid, it 
experiences the Magus force. This force has the same ef- 
fect on the vortex as magnetic field on a moving charged 
particle, so that it can be implemented through an effec- 
tive gauge field ^ in the vortex Hamiltonian: 



(39) 



2mv 2 

where = 1), for a Galilean-invariant superfluid, 

V X ^ = 27rpsi; , (40) 

and ps is the superfluid density (density of Cooper pairs) . 
On a lattice, recent work^ has argued that the rhs has to 
be replaced by difference between the density of Cooper 
pairs in the superfluid and half the density of electrons 
in a proximate solid. The presence of the Magnus force 
defines a corresponding cyclotron frequency: 



27rps 



(41) 



We assume that the vortex trap is isotropic, defining a 
single trap harmonic frequency lo^. Hence, both the trap 
and the Magnus force support circular classical trajec- 
tories of the vortex, making it easy to diagonalize the 
Hamiltonian H39I) : 



(42) 




e/cov 



e/cov 



e/cOv 



FIG. 7: (color online) The full LDOS p = po + Pi at the 
vortex center as a function of energy and quasiparticle cut-off 
A. The dotted red curve is po(e), which excludes the resonant 
scattering and does not depend on any finite cut-off. 



with 



Wo 



(43) 



Wo 



L0+ = Wo + — 



and 



V 2 

TOvWq 



Px 



V2 



^ Pxv 



V2 



(44) 



iPyv 



V2 



rrivLUo 



V2 



The full Hamiltonian includes coupling between the vor- 
tex and quasiparticles: 



H = H. 



E 

nodes 



(45) 



and we can proceed solving it just like before. 

The Magnus force introduces only a few small changes 
in the LDOS obtained so far. First, the zero-point quan- 
tum motion of the vortex is now controlled by the fre- 
quency Wo, instead of Wv, which can be seen from 144|l . 
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FIG. 8: The full LDOS at the vortex center as a function of 
energy, for different values of the Magnus cyclotron frequency 
i^m (^Ea. I41|l . The plots are offset vertically for clarity. 



The vortex ground-state wavefunction is still isotropic, 
but extends in space to distances ~ (mv^o)"^^^ from 
the trap center. Second, there are two small parameters, 
a+ — y/m^ujQ I a;+ and a_ — ^Jrru^UJ^ j lo - that shape the 
perturbation theory. The frequencies a;+ and w_ are im- 
portant for the resonant scattering of quasiparticles from 
the vortex, and in principle, discontinuities in the LDOS 
can occur at these frequencies. The Magnus force lifts 
the degeneracy of vortex eigenmodes in an isotropic trap, 
so that the lowest excited states correspond to left and 
right handed circular motion of the vortex, with different 
energies. 

In the Fig. |Hlwe compare the quasiparticle LDOS for 
several values of the Magnus cyclotron frequency ciJni- In 
general, no qualitative changes arise when is finite, 
provided that it does not become too large [uj^ ^ 2ujy <^ 
a- would invalidate the perturbation theory). As 
(jJin grows, the sub-gap peak shifts toward smaller ener- 
gies. Somewhat surprisingly the sharp features, such as 
the LDOS discontinuities, seem to become smoother as 
Wni grows, although they do not completely vanish. In 
general, there is only one discontinuity at e = 



IV. DISCUSSION AND CONCLUSIONS 

We have explored how the quantum motion of vortices 
affects quasiparticle spectra in clean d-wave supercon- 
ductors at zero temperature. Keeping only the statisti- 
cal interaction between vortices and quasiparticles in a 
continuum Bogoliubov-de Gennes model, we have found 
that quantum oscillations of a vortex could be responsible 
for observation of the "core states" in the STM experi- 
ments on cuprates. The emerging physical picture is that 
there are no bound states in the vortex core of an ideal d- 
wave superconductor, but quasiparticles can instead ex- 
perience resonant scattering from a vortex. Vortices are 
localized either by their neighbors in a vortex lattice, or 
by pinning impurities, and execute zero-point quantum 
oscillations in their harmonic traps. When a quasiparti- 
cle scatters from a vortex, it can excite a virtual higher- 
energy state of vortex oscillations, and effectively form a 
short-lived metastable bound state with the vortex. This 
leads to a peak of the local density of states at the appro- 
priate resonant energy in the vicinity of the vortex core. 
The peak lies inside the superconducting gap because the 
vortex zero-point motion extends over length-scales that 
are significantly larger than the coherence length. 

Our numerically calculated LDOS has many similar- 
ities with experimental observations at energies smaller 
than the superconducting gap. The low energy scans 
of LDOS at various distances from the vortex core look 
qualitatively the same as the experimental measure- 
ments. We do not obtain a zero-energy peak that was 
originally predicted by meanfield BCS calculations^, and 
we find a small sub-gap peak that gradually vanishes 
with growing distance from the vortex core. We also 
find other features in the LDOS energy scans, such as 
discontinuities and secondary peaks and dips, but these 
may be blurred and too small to observe in realistic cir- 
cumstances. According to our model, energy of the sub- 
gap peak turns out to be an increasing, but not a lin- 
ear function of the superconducting gap. Experimental 
data appear to be consistent with the linear dependence, 
but due to large error margins in the measurements, and 
many simplifications built in our model, we believe that 
this detail cannot rule out the great importance of vor- 
tex quantum motion for quasiparticle spectra in d-wave 
superconductors. 

The ability to reach a comprehensive physical picture 
of the sub-gap quasiparticle spectra comes with a price: 
our model contains many simplifications of realistic cir- 
cumstances. Most approximations are controlled in the 
sense that there is a limit in which they become accurate. 
Away from that limit there will be quantitative modifica- 
tions of the results due to other effects, such as thermal 
fluctuations, decoherence, and possible competing orders 
inside the vortex cores. We have also ignored anisotropies 
of d-wave superconductors. Perhaps the least justified 
simplification is inclusion of only the statistical interac- 
tions between vortices and quasiparticles. The missing 
Doppler shift presents an effective potential to quasipar- 
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tides near the vortex core, and may enhance the resonant 
scattering. We do not expect that it could quaUtatively 
change the results, but its quantitative contribution need 
not be small. 

We also reiterate here the original motivation for con- 
sidering the vortex zero-point motion model. This arose 
from its consequences in the spatial dependence of the 
LDOS: it was arguedSii that the Aharanov-Bohm phases 
acquired by vortices from the background density of 
electrons lead naturally to periodic modulations in the 
LDOS. This was proposed as an explanation for the mod- 
ulations observed in STM experiments^. Here we have 
shown that the same model can also help account for 
some crucial aspects of the energy dependence of the 
LDOS. Further, as argued in Ref. 3, there is also a quan- 
titative consistency between the energy scales deduced 
from the spatial and energy dependencies. From the ex- 
tent of the spatial modulationsj^ and an estimate of the 
trapping potential on a vortex, we were able to obtain an 
estimate of the vortex oscillation frequency.'^ This num- 
ber is consistent with the observed positio n'^i^^i^^ of the 
sub-gap peaks in the energy dependence of the LDOS. 

Future experiments may provide more ways to detect 
consequences of vortex quantum motion in d-wave super- 
conductors. Neutron or light scattering measurements 
could provide a direct observation of the vortex oscilla- 
tions. The spatial extent of the LDOS modulations pro- 
vides means for measuring one important energy scale, 

hajviTivVf , where Wv is the harmonic frequency of a lo- 
calized vortex. Our present prediction is that the sub- 
gap peak in the quasiparticle LDOS occurs at the en- 
ergy proportional to this scale. An increased sensitiv- 
ity of STM measurements in the future could allow re- 
solving the jumps in the LDOS at energies that reflect 
the discrete spectrum of the trapped quantum vortex. 
This would provide means for measuring another impor- 
tant energy scale, Combined knowledge of both en- 
ergy scales reveals the vortex mass mv and the strength 
of the vortex trapping force (inter- vortex interactions), 
which in turn can be related to other parameters, such 
as the superconducting gap and magnetic field. Effects 
of the Magnus force reduce to quantitative modifications 
of these energy scales, as discussed in the Section UlI Dl 
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oscillator into (|HJ , and have the same form as the Hamil- 
tonian (O in which the vector operator a is substituted 
by the vector c-number: 



Q'(nix,ni„),(n2x,n2y) (^) 



27r 



z X [r - Ty) _ 



y/nixln2xlniyln2yl 

(Al) 



" 2|r-rv|2 

Hn{x) are Hermite polynomials of order n. One way to 
analytically calculate these integrals is to use the follow- 
ing formula: 



^^i-^ = -.x(Vlog|r-r.|) 

and expand: 

log \r - Tvl = log \z - Zvl = log 
1 



(A2) 



1-^ 



log|z>| 



■log 1 



■ log 1 



= log|z>|-Re5:i(^) 



(A3) 



Coordinates have been represented by complex numbers 
z and Zy] z> and z^ are the greater and lesser by modulus 
of z and z^ respectively. The effective gauge-field (|A1|I 
is: 



<l(r!.ix,nia),("2x,n2H) (''") ^ 



(A4) 



(X 2 X V / d Zv 



log|z; 



n=l \ ^ 



— nivt^vl^vl^ rr TJ M TJ 

The Hermite polynomials are easily written in terms of 
Zv and z* at low orders n. For any particular set of 
{nlx^ niy, n2x, n2y) only one term in the expansion above 
gives a non-zero contribution upon integration of the 
phase of z^. For example, when nix = niy = n2x = 
= all Hermite polynomials are equal to unity, and 
the only term that does not depend on phase is the log- 
arithm: 



1(0,0), (0,0) 



i; X V 



d|zv|e-'"-'^-'^-l'log(|z|) + 



APPENDIX A: CALCULATION OF VERTEX 
OPERATORS 



Here we derive general expressions for the "matrix ele- 
ments" (jSJ, which play the role of vertex operators in the 
perturbation theory. The matrix elements are evaluated 
by substituting the wavefunctions of the linear harmonic 



dlzvle-^'^-l^-I'lo: 



(A5) 



m^jUJ^z X r—\ogir) / dr^r^e 
dr Jq 



zx — [I- e-"-"-'- 
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This is the result written in A similar calculation 

produces the expression (j^ . 

APPENDIX B: CALCULATION OF THE 
ONE-LOOP CORRECTION TO THE LDOS 

Here we substitute (|28ll and (|27|) into H35|l and derive 
the lowest order correction to the LDOS. We systemati- 



cally use the formula: 

1 P 

-— = ^7TS{u;-e) (Bl) 

w — e + tv^ UJ — e 

to extract the imaginary part in H35|l . After a few alge- 
braic manipulations one obtains the following correction 
to the LDOS: 



h,l2 I I 

-TT^T (r)f^^^'^ /fie(6-^.) 130 \ , fm{e) fiO \ t 

""^hM^^'y fiO Q-e{-e) J ^h,\e\;L\e\~u;,.y QQ Q ^ - Lu^) J ^iM-'^-^'i2M ^ fiO & - e{-e) J ^2,^1^^ > 

+ y dfcidfc2 J ;^ ,^(rj fiO Q-^J^hM;iM-^.[ fiO Q - ei-e - u;^) ) ^'M-^^'^hM [ fiO fi-^ j ^'2,^2^^^ 

+ /dfc,dfcT Ml'^^^^) V- (r) 

+ j dfe2dfcJ,^l^|(r) fiO Q-&{-e) fio fi^^^ j ^''^^'^^^^ 1^ ^0 fi^ j ^'^.fc.i^) 

+ y dfcdfci J ;^,,^irj fiO j ^'i.fei;',^- fiO Q-^^ )^iMhM[ fiO fi - e(-e) j ^'2,|e|l^J j • 



This expression for pi is written assuming that the 
sample is infinite. In order to carry out numerical cal- 
culations it is necessary to impose a finite sample ra- 
dius R and thus quantize energy levels of the unper- 
turbed Hamiltonian. The quantized energies take values 
e„ = ±fc„, where A;„i? -|- 5 is the n*'' zero of the Bessel 
function {S is a phase shift). In the region ^ 1, 
the energy levels are separated hy n/R, and this allows a 
simple rule of thumb for conversion between the continu- 
ous and discrete computations. First, the infinite-sample 
wavefunctions need to be renormalized to unity on the fi- 
nite sample. Up to a small error of the order (fc„i?)~^ this 
amounts to multiplying the wavefunction by y^ir/R. Sec- 
ond, all integrations over continuous radial wave-vectors 
k in the expression HB2() must be replaced by summations 
over discrete fc„, with a measure factor of n/R. Finally, 
the principal part P/x is taken to deviate from 1/x in 
the interval of width tt/R from x = 0, passing smoothly 
through zero at a; = 0. 

As a result of energy discretization, the LDOS is 
strictly speaking defined only at discrete energy levels 
e„. Such discretization is naturally lifted if there are 
additional fluctuations in the problem that broaden the 



energy levels beyond ir/R (for example, due to finite 
temperature or disorder). However, the level-broadening 
mechanisms are not important in the R 00 limit, 
which is reached in numerics if the majority of states 
satisfy fc„i? ^ 1 (that is, if the energy cut-off A is large, 
AR ^ 1). The typical values that we used in numerics 
were R = 200, A = 4, Wv = 1, and -5<l<5. 



APPENDIX C: MICROSCOPIC THEORY FOR A 
d-WAVE SUPERCONDUCTOR 

This appendix connects the results of the present paper 
to those of paper I. Here we will extend the numerical 
results of I to d-wave superconductors. However, because 
of the gradient expansion necessary in this method, we 
are unable to extend these results to very small vortex 
core sizes. Consequently, we will not find an elimination 
of the zero bias peak in the LDOS. 

As discussed in I, we will use the gap operator given in 
Eq. (3) of I, with A(r) ~ A(r)e'^. In polar coordinates 
which are more convenient for our problem here this gap 
operator reads 



13 



A =^T2-T [(3A - 3rA' + r^A") + 4(-2rA + r^A')dr + 4i(-3A + rA')de + AAir^d^ - + 2irdr 



(CI) 



As a generalization of Eqs. (21) and (31) of I we write 
the quasi-particle wave functions as 




where the 

(j>,n,j{r) = (l)m,j{r, 9) = exp[-im6']0„j(r)/V27r (C3) 

[which generahze Eq. (29) of I] now form a complete set of 
eigenfunctions to the kinetic energy operator and satisfy 
the normalization condition 

(C4) 

JQ JQ 

I 



The Bogoliubov-de Gennes equations reduce to a matrix 
equation quite similar to Eq. (32) of I. The essential dif- 
ference is that there is now an infinite hierarchy of cou- 
pled angular momentum channels and the matrices 
and A now have matrix elements 

'^mj-m'j' = T (^"^f^ ~ ^m.m'Sjj' , (C5) 

/•2-lT I'Rq 

Amj-m'j' ^ d9 drr(f>l^ j{r,e)A(t)m,_j,{r,9). 
Jo Jo 

(C6) 

Using the gap operator given in Eq. IjCip we can explicitly 
evaluate A„y .„j/j' and obtain 



A™,;™-,v = S„,,^,n+3 f " dr [(-4(a™+3jvr/i?o)2 + 3) A - 3rA' + r^A"] cbm+3.j'{r) 

Jo oKpr 

-f4[(2m + 3)A + rA'] 0-+3,,,(r)} 
+ <5,„',™-i / " dr ^^^{ [(-4(a™_i,,vr/i?o)' - 1)A + rA' + r^A"] ^m-i.^r) 

Jo oKpT 

+ 4[(2m-l)A-rA']0+_,_^.,(r)} . (C7) 



The functions 

= ^"'"'J'' Jm±i{am.jr/RQ) (C8) 
n.Q\Jm+iyam.j)\ 

are modifications of the normalized Bessel functions 
(pm .j{r) defined in Eq. (29) of I which arise when tak- 
ing derivatives of (ftrnjir)- According to Eq. (jC7|) to' has 
to equal to — 1 or to -|- 3 for a matrix element not to vanish. 
The reason for this is that while the vortex in the order 
parameter mediates a change in to (minus angular mo- 
mentum) by 1 when scattering the hole-like component 
of the quasi-particlc into an electron- like component, due 
to the d-wave symmetry of the gap operator we have an 



additional change by ±2. As a consequence of this, the 
eigenvalue problem breaks up into four apparently in- 
dependent eigenvalue problems which we all truncate at 
large angular momentum m,nax- 

Let us say that the eigenvector Vf^ lies in the sector 
s{£) = 0, 1, 2, or 3 if '^i contains only components c^u+s 
and (i4^+s_i with v an integer. It is possible to write the 
eigenvalue problem for each of these four sectors in terms 
of a real and symmetric band-diagonal matrix. For con- 
creteness, let us consider the sector s — Q with TO^ax — 5 
(which, of course, for practical calculations should be cho- 
sen much larger). We then have 
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Vo 



A_4;-: 

Aq.-I 



A^: 
AL 



\ 



Ao;3 

T+ 

A4.3 



(C9) 



Here, and Am;m' are themselves TVq x iVg matrices 
whose matrix elements are given above. We can now 
relate each eigenstate in sector 1 to one in sector and 
each eigenstate in sector 3 to one in sector 2 by making 
use of the property that if [^^(r), (r)]^ is an eigenstate 
of HBdG 'with eigenvalue e^, then [— u|(r), u|(r)]-^ is an 
eigenstate of Hg^jQ with eigenvalue —ei. This property 
implies the symmetry 
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which we have checked explicitlyii^ In terms of the 
eigenvectors ^ = (c4iy+s7 '^4i/+s-i) we find that if ^' 
is an eigenvector in sector s with eigenvalue e, then 
* = (c4^+(i_s),(i4i.+(i-s)-i) with Q 



d-m,j and 



is an eigenvector with eigenvalue — e in 
sector 1 — s (mod 4) . It is therefore sufficient to only con- 
sider the sectors and 2 and make use of this property 
to obtain all eigenstates. 

An interesting property of having no mixing between 
the four sectors is that the quasi particle amplitudes 



^[''•'(r) and vf{r) can be written as 
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Physically measurable quantities like the local DOS 
which only involve |M^''^(r)p (and |z;^*''(r)p) are there- 
fore invariant under all symmetry operations expected 
for a (ia,2_j,2-wave order parameter, i.e. rotations by 90° 
as well as reflections along the nodal or anti-nodal di- 
rections. Before considering the effect of the zero-point 
fluctuations of the vortex on the local DOS let us first 
have a look at the case of a static vortex which was al- 
ready studied in Refs. 0|3 One advantage of using the 
gap operator given in Eq. (|C1|) is that most of its matrix 
elements Amj;m'j' cancel. We can therefore calculate the 
local DOS for reasonably large system sizes with moder- 
ate computational effort. 



In Fig. 1^1 we show results for the angular average of 
the tunneling conductance G{r,uj) for various values of 
r calculated for a system with radius kpR = 120 and a 
cutoff rriniax = 60 for the angular momentum channels. 
Directly at the vortex center there is a peak in the local 
DOS near the Fermi level which is due to a continuum of 
unlocalized states whose envelopes take there maximum 
near the vortex center but which have tails leaking out 
in the nodal directions. These results are consistent with 
those found in Refs.|6ii2i. As we move away from the vor- 
tex center the peak in the local DOS gradually vanishes 
and the tunneling conductance assumes the familiar fea- 
tures of a d-wave bulk spectrum (at finite temperature). 

To calculate the transition matrix elements M"^, which 
were defined in Eq. (15) of I we will use the Hellmann- 
Feynman theorem which in the present case of a d-wave 
gap operator reads 
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The first two terms are as in the s-wave case but since 
our d-wave operator involves two derivatives, partial in- 
tegration now also leads to a boundary term involving 
the bulk gap. In the s-wave case angular momentum 
was a good quantum number and for a transition ma- 
trix element not to vanish we had the selection rule that 
the angular momenta of the two wave functions had to 
differ by one. As already emphasized above, in the d- 
wave case angular momentum is not a good quantum 
number. Instead, each eigenstate belongs to one of four 
sectors which only contains certain angular momentum 
channels. As a result of this we find that the transi- 
tion matrix elements M^^, = ^ are only non-zero if 
s'{£') — s{£) + 1 (mod 4). In this case we can again ex- 
press all wave functions in terms of their Fourier-Bessel 
components and obtain 
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FIG. 9: Angular average of the tunneling conductance G — 
dl /dV for a superconductor with d^2_y2 wave symmetry and 
a static vortex at the origin as a function of ui. We have 
chosen kpS, ~ 2.5 and T = 0.02 Ep. The upper curve is for 
r = and the curves below are for kpr = 4, 8, . . . , 20. 



FIG. 10: Angular average of the tunneling conductance G = 
81/ dV for a superconductor with d^2_y2 wave symmetry and 
a vortex with vortex mass m„ = as a function of oj. We 
have chosen fc^C = 2.5 and T = 0.02 i?_F. The upper curve is 
for r = and the curves below are for kpr = 4, 8, . . . , 20. 



where the matrix elements /C^y and 'C^™^ were already 
defined in Eqs. (40) and (41) of I and rin{m) = 1 — 25n,m 
leads to an extra minus sign if m ~ n. 

The tunneling conductance calculated for the param- 
eters used above but with a vortex mass equal to the 
mass of an electron is shown in Fig. 1101 As seen in experi- 
ments, the central peak in the vortex center is suppressed 
and weight is shifted to both sides away from the Fermi 
level. Although the mechanism for redistributing spec- 
tral weight is the same as in the s-wave case no satellite 
peaks are visible in the d-wave case. The redistribution 
of spectral weight is further illustrated in Fig. 1111 where 
we compare our results for a static vortex with those for 



a moving vortex. As can be seen in that figure, choosing 
a smaller vortex frequency ujy leads to a much stronger 
redistribution of spectral weight. This is a direct conse- 
quence of Eqs. (16) and (17) of I. Also, although there are 
no sub-gap peaks visible in the tunneling conductance, 
the curves for the difference between the tunneling con- 
ductance for a moving and that for a static vortex show 
peaks at an energy significantly smaller than the gap en- 
ergy, which according to Eq. (16) of I should be of the 
order of the vortex frequency. 

Although our theory explains the redistribution of 
spectral weight from the zero bias level towards higher 
energies, it cannot explain the exact line shape as ob- 
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FIG. 11: Comparison of the tunneling conductance G — 
dl/dV for a superconductor with d^2_y2 wave symmetry in 
the presence of a static or moving vortex. The upper curve 
shows the tunnehng conductance at the vortex center for a 
static vortex (with kp^ = 2.5 and T — 0.02 Ep). As we allow 
for a finite vortex mass and vortex frequency the cental peak 
decreases and weight is shifted away from the Fermi level. 
This is shown in the two curves below for which we have cho- 
sen m„ = me and tj„ — 0.2 Aq or = 0.5 Ao respectively. 
Finally, the two curves at the bottom of the figure represent 
the difference between the results for the case of a moving 
vortex and that of a static vortex. 



an unphysical peak at « +Ao as seen for rriyaax = 24 or 
36 is basically absent for TOmaa; = 60.) 

The above scenario is quite different from the case of 
a static vortex where a cutoff of rritnax — 24 suffices to 
obtain convergence. We therefore conclude that it must 
be the transition matrix elements M^^, which are very 
sensitive to the boundary conditions imposed by the cut- 
off mraax ■ Another restriction on the applicability of our 
theory is our expansion of the gap operator A(r — R(t)) 
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FIG. 12: Comparison of results for the tunneling conduc- 
tance G = 81/ dV at the vortex center of a superconductor 



served in experiment. First of all, the system sizes stud- 
ied are much smaller than in the s-wave case, leading to 
a poorer spectral resolution. Angular momentum is not 
a good quantum number and a large number of angu- 
lar momentum channels is needed to obtain convergence. 
This is shown in Fig. 1121 where we compare results for 
the tunneling conductance at the vortex center, calcu- 
lated with different cutoffs m^ax implemented. As can 
be seen in that figure, we have not even achieved conver- 
gence for rTimax = 60, the value used in all the previous 
figures. However, it seems that we are not too far away 
from the m^ax ~^ oo limit. (It should also be noted that 



tex. As in the previous figures, we have considered a sys- 
tem with radius kpRo = 120, kp^ = 2.5 and T = 0.02 Sf. 
The only difference is that we have used different cutoffs 
(jTimax ~ 24, 36, 48, and 60) for the angular momentum chan- 
nels. 



to leading order in R(r). Since the quantum zero-point 
motion of the vortex basically extends roughly over a 
distance of l/y/rrivUJv from the vortex center this ap- 
proximation is only expected to be reasonable as long 
as ^ 1/y^to^A^. We are therefore not able to study 
systems with very small coherence lengths. 
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